Solving Stochastic Differential Equations with Jump-Diffusion 
Efficiently: Applications to FPT Problems in Credit Risk 

Di Zhang 1 and Roderick V.N. Melnik 1 

1 Mathematical Modelling and Computational Sciences, Wilfrid Laurier University, 
Waterloo, ON, Canada N2L 3C5, Email: rmelnik@wlu.ca 



Abstract. The first passage time (FPT) problem is ubiquitous in 
many applications. In finance, we often have to deal with stochastic 
processes with jump-diffusion, so that the FTP problem is reducible 
to a stochastic differential equation with jump-diffusion. While the 
application of the conventional Monte-Carlo procedure is possible 
for the solution of the resulting model, it becomes computationally 
inefficient which severely restricts its applicability in many practi- 
cally interesting cases. In this contribution, we focus on the develop- 
ment of efficient Monte-Carlo-based computational procedures for 
solving the FPT problem under the multivariate (and correlated) 
jump-diffusion processes. We also discuss the implementation of 
the developed Monte-Carlo-based technique for multivariate jump- 
diffusion processes driving by several compound Poisson shocks. 
Finally, we demonstrate the application of the developed method- 
ologies for analyzing the default rates and default correlations of 
differently rated firms via historical data. 

Keywords. Default Correlation, First Passage Time, Multivari- 
ate Jump-Diffusion Processes, Monte-Carlo Simulation, Multivari- 
ate Uniform Sampling Method. 

AMS (MOS) subject classification: 60H35, 65C05, 68U20. 

1 Introduction 

Many problems in finance require the information on the 
first passage time (FPT) of a stochastic process. Mathe- 
matically, such problems are often reduced to the evalu- 
ation of the probability density of the time for a process 
to cross a certain level. Recent research in finance the- 
ory has renewed the interest in jump-diffusion processes 
(JDP), and the FPT problem for such processes is appli- 
cable to several finance problems, such as pricing barrier 
options [TJ [5] , credit risk analysis [3] . While in other ar- 
eas of applications the FPT problem can often be solved 
analytically, in finance we usually have to resort to the 
application of numerical procedures, in particular when 
we deal with jump-diffusion processes. 

Among numerical procedures, Monte-Carlo methods 
remain a primary candidate for applications. However, 
the conventional Monte-Carlo procedure becomes com- 
putationally inefficient when it is applied to the jump- 
diffusion processes. Many researchers have contributed to 
the field of enhancement of the efficiency of Monte-Carlo 
simulations. Atiya and Metwally [HE] have recently de- 
veloped a fast Monte-Carlo-type numerical method to 
solve the FPT problem in the one-dimensional case. 

Note that apart from the pricing and hedging of op- 
tions on a single asset, practically all financial applica- 



tions require a multivariate model with dependence be- 
tween different assets. In particularly, jumps in the price 
process must be taken into account in most of the appli- 
cations [5l [6]. In this contribution, we generalize our pre- 
vious fast Monte-Carlo method (for non-correlated jump- 
diffusion cases) to multivariate (and correlated) jump- 
diffusion processes. The developed technique provides an 
efficient tool for a number of applications, including credit 
risk and option pricing. We also discuss the implementa- 
tion of the developed Monte-Carlo-based techniques for 
a subclass of multidimensional Levy processes with sev- 
eral compound Poisson shocks. Finally, we demonstrate 
the applicability of this technique to the analysis of the 
default rates and default correlations of two different cor- 
related firms via a set of empirical data. 

The article is organized as follows: Section [5] describes 
the mathematical model. The methodologies are pre- 
sented in Section [31 applications and discussions are given 
in Section [5J Concluding remarks are given in Section [5J 

2 Mathematical model 

In this section, first we present a probabilistic description 
of default events and default correlations in finance. Next, 
we describe the multivariate jump-diffusion processes and 
provide details on the first passage time distribution un- 
der the one-dimensional Brownian bridge (the algorithms 
which is used to generate correlated first passage times 
will be described in Section [3]). Finally, we present kernel 
estimation in the context of our problem that can be used 
to represent the first passage time density function. 

In finance, a firm i defaults when it can not meet its 
financial obligations, or in other words, when the firm 
assets value Vi (t) falls below a threshold level Dy t it) . In 
this contribution, we use an exponential form defining 
the threshold level Dy^t) — Ki exp(jit) as proposed by 
Black and Cox [7], where 74 can be interpreted as the 
growth rate of firm's liabilities. Coefficient m captures 
the liability structure of the firm and is usually defined as 
a firm's short-term liability plus 50% of the firm's long- 
term liability [8]. If we set Xi(t) = ln[Vj(i)], then the 
threshold of Xi(t) is Di(t) = jit + ln(/«i). Our main 
interest is in the process Xi(t). 

In the market economy, individual companies are in- 
evitably linked together via dynamically changing eco- 
nomic conditions [5]. Therefore, the default events of 
companies are often correlated, especially in the same in- 



dustry. Take two firms i and j as an example, whose stochastic differential equation: 
probabilities of default are Pi and Pj , respectively. Then 
the default correlation can be defined as 
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PiPi 



Pij 



(i) 



where P« is the probability of joint default. 

Zhou [5] and Hull et al [9] were the first to incorpo- 
rate default correlation into the Black-Cox first passage 
structural model. They have obtained the closed form so- 
lutions for the joint probability of firm 1 to default before 
T\ and firm 2 to default before T2. However, none of the 
above known models includes jumps in the processes. At 
the same time, it is well-known that jumps are a major 
factor in the credit risk analysis. With jumps included in 
such analysis, a firm can default instantaneously because 
of a sudden drop in its value which is impossible under 
a diffusion process [3]. Therefore, for multiple processes, 
considering the simultaneous jumps can be a better way 
to estimate the correlated default rates. 

A natural approach to introduce jumps into a multi- 
dimensional model is to utilize the compound Poisson 
shocks. The dates of market crashes can be modeled as ar- 
rival times of a standard Poisson process Nt , which leads 
us to the following model for the log-price processes of d 
assets [5]: 



dXi — [Lidt + <TijdWj + dZi 
j'=i 

= mdt + a-idWi + dZ h (4) 
where Wi is a standard Brownian motion and cr, = 

We assume that in the interval [0, T], the total num- 
ber of jumps for firm i is Mj. Let the jump instants be 
T\,T%, ■ ■ ■ , 7V 4 . Let To = and Tm;+i = T. The quan- 
tities Tj equal to interjump times, which are Tj — Tj-\. 
Following the notation of |4J, let Xi(T~) be the process 
value immediately before the jth jump, and Xi(T+) be 
the process value immediately after the jth jump. The 
jump-size is Xi(Tj) — Xi{T~) and we can use such jump- 
sizes to generate Xi(Tf) sequentially. 

Although for jump-diffusion processes, the closed form 
solutions are usually unavailable, yet between each 
two jumps the process is a Brownian bridge for one- 
dimensional jump-diffusion process. Let B(s) be a Brow- 
nian bridge in the interval [Tj-i,Tj] with B{T+_ X ) = 
XiiTf^) and B{Tr) = X t (Tp. If JQ(Tr) > £><(<), 
then the probability that the minimum of B(si) is always 
above the boundary level is [4] 



Xi(t) = mt + Bt(t) + Zt(t), i = 1,2, • • • , d, 



where B(t) is a d-dimensional Brownian motion with co- 
variance matrix a = (<7y) which can be written as 



2{X i (TT 1 )-D i (t)][X i (T )-D i {t)] 



B i {t) = Y,°v W M 



and Wj(t) is the standard Brownian motion. For i-th 
asset, {Yij}J^ 1 are i.i.d. d— dimensional random vectors 
which determine the sizes of jumps in individual assets 
during a market crash. At the j-th shock, the jump-sizes 
of different assets Yij may be correlated. 

This model contains only one driving Poisson shock 
which stands for that the global market crash affecting 
all assets. Sometimes it is necessary to have several inde- 
pendent shocks to account for events that affect individ- 
ual companies or individual sectors rather than the entire 
market. In this case we need to introduce several driving 
Poisson processes into the model, which now takes the 
following form [5]: 

m w t * 

Xi(t)=Hit + Bi(t)+J2J2 Y ijk, * = 1, 2, • • • ,d, (3) 
k=ij=i 
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(5) 



This implies that B(si) is below the threshold level, which 
means the default happens or already happened, and its 
probability is 1 — Pij. 

For firm i, after generating a series of first passage times 
Si, we use a kernel density estimator with Gaussian kernel 
to estimate the first passage time density (FPTD) /. The 
kernel density estimator is based on centering a kernel 
function of a bandwidth as follows: 



1 N 

— YK(h,t 



(G) 
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The optimal bandwidth in the kernel function K can be 
calculated as HOI: 



h 



Op I 



= [2N^ [ 
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{ft'fdt 



(7) 



where N is the number of generated points and ft is the 
true density. Here we use the approximation for the dis- 
tribution as a gamma distribution as proposed in [4]: 



ft 



^ _1 exp(-ai). 



(8) 



where , • • • , N™ are Poisson processes driving m inde- 
pendent shocks and Yy/c is the size of jump in i-th com- 
ponent after j-th shock of type k. The vectors {Yijk}f=\ 
for different j and/or k are independent. 

First, let us consider a firm i, as described by Eq. 



3 Methodology of the solution 



Let us recall the conventional Monte-Carlo procedure in 
application to the analysis of the evolution of firm Xi 
within the time horizon [0, T]. We divide the time hori- 
@, such that its state vector Xi satisfies the following zon into n small intervals [0, fi], [^1,^2], [t n ~i,T] as 



shown in Fig. [ija). In each Monte-Carlo run, we need 
to calculate the value of Xi at each discretized time t 
without explicitly distinguishing the effects of the jump 
and diffusion terms [11] . As usual, in order to reduce dis- 
cretization bias, the number n must be large 12J. The 

Conventional Monte- Carlo method 
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Figure 1: Schematic diagram of (a) the conventional 
Monte-Carlo and (b) the uniform sampling (UNIF) 
method. 

conventional Monte-Carlo procedure exhibits substantial 
computational difficulties when applied to jump-diffusion 
processes. Indeed, for a typical jump-diffusion process, 
as shown in Fig. QJb), let 3j-i and Tj be any successive 
jump instants, as described above. Then, in the conven- 
tional Monte-Carlo method, although there is no jump 
occurring in the interval \Tj-\,Tj], yet we need to evalu- 
ate Xi at each discretized time t in [Tj_x, Tj]. This very 
time-consuming procedure results in a serious shortcom- 
ing of the conventional Monte-Carlo methodology. 

To remedy the situation, two modifications of the con- 
ventional procedure were recently proposed [TJ [1] that al- 
low us a potential speed-up of the conventional method- 
ology in 10-30 times. One of the modifications, the uni- 
form sampling method (UNIF), involves sampling using 
uniform distribution. The other, inverse Gaussian density 
sampling is based on the inverse Gaussian density method 
for sampling. Both methodologies were developed for the 
univariate case. 

The major improvement of the uniform sampling 
method is based on the fact that it only evaluates Xj at 
generated jump times, while between each two jumps the 
process is a Brownian bridge (see Fig. [11(b)). Hence, we 
just consider the probability of Xi crossing the threshold 
in Tj) instead of evaluating Xi at each discretized 

time t. More precisely, in the uniform sampling method, 
we assume that the values of X^X^-J and Xi(Tj~) are 
known as two end points of the Brownian bridge, the 
probability of firm i defaults in (Tj_i, Tj) is 1 — which 
can be computed according to Eq. ([5]) . Then we generate 
a variable Sj from a distribution uniform in an interval 
\Tj—i,Tj—x + Jp.T 1 ]■ If the generated point Si falls in 
the interjump interval [Tj—i, Tj], then we have success- 
fully generated a first passage time Si and can neglect the 
other intervals and perform another Monte-Carlo run. On 
the other hand, if the generated point Si falls outside the 
interval \Tj-\,Tj] (which happens with probability P%j), 
then that point is "rejected" . This means no boundary 
crossing has occurred in the interval, and we proceed to 
the next interval and repeat the whole process again. 



In what follows, we focus on the further development 
of the uniform sampling method and extend it to mul- 
tivariate and correlated jump-diffusion processes. First, 
we consider there is only one driving Poisson shock with 
arrival rate A as described in Eq. ([2]). The distribution 
of (Tj — Tj_i) are the same for each firm. The jump- 
size can be generated by a given distribution which can 
be different for different firms to reflect specifics of the 
jump process for each firm. In order to implement the 
UNIF method for the multivariate processes, we need to 
consider several points: 

1. We exemplify our description by considering an expo- 
nential distribution (mean value /it) for (Tj — Tj_i) 
and a normal distribution (mean value \xj and stan- 
dard deviation aj) for the jump-size. We can use 
any other distribution when appropriate. 

2. An array IsDefault (whose size is the number of 
firms denoted by -/Vfi rm ) is used to indicate whether 
firm i has defaulted in this Monte-Carlo run. If the 
firm defaults, then we set IsDef ault(i) = 1, and will 
not evaluate it during this Monte-Carlo run. 

3. Most importantly, as we have mentioned before, the 
default events of firm i are inevitably correlated with 
other firms, for example firm i + 1. Hence, firm 
j's first passage time Si is indeed correlated with 
Si+i - the first passage time of firm i + 1. Wc 
must generate several correlated Si in each interval 
[Tj_i, Tj-i + T \ T p7. x ] which is the key point for mul- 
tivariate correlated processes. 

Note that each process is a Brownian motion in the 
interval [Tj-i, Tj], so we can compute the correlation 
coefficient pis+i of firms i and i + 1 by using Zhou's 
model without jumps [5] and then use this value for 
modeling correlated Si and Sj+i. Let us introduce 
a new variable bij = J' 1 . Then we have Si = 
bijYi + Tj_x, where are uniformly distributed in 
[0, 1]. Moreover, the correlation of Yi and Y i+1 is the 
same as Si and s^+i. The correlated uniform random 
variables Y\, Yi, ■ ■ ■ can be generated by using the 
sum-of-uniforms (SOU) method |13j . 

Next, we will describe our algorithm for multivariate 
jump-diffusion processes, which is an extension of the one- 
dimensional case developed earlier by other authors (e.g. 

[DEI). 

Consider iVg rm firms in the given time horizon [0, T]. 
First, we generate the jump instant Tj by generating in- 
terjump times (Tj — Tj_i) and set all the IsDef ault(i) = 
to indicate that no firm defaults at first. 

From Fig. [TJb) and Eq. (TJ|, we can conclude that for 
each process Xi we can make the following observations: 

1. If no jump occurs, as described by Eq. ([4]), the in- 
terjump size (Xi(T^) — Xi(Tj_-y)) follows a normal 
distribution of mean yn(Tj — Tj-i) and standard de- 
viation Oi^jTj — Tj-i. We get 

Xi(Tf) ~ X t (T+_ x ) + Wj + CjAT (0, Tj), 

where the initial state is X(0) = Xi(T^). 



2. If jump occurs, we simulate the jump-size by a nor- 
mal distribution or another distribution when appro- 
priate, and compute the postjump value: 

X l (T+) = X l {T-) + Z l {T j ). 

This completes the procedure for generating before- 
jump and postjump values Xi(Tj~) and Xi(Tj ). As be- 
fore, j = 1, • • • , M where M is the total number of jumps 
for all the firms. We compute i^y according to Eq. ([5]). 
To recur the first passage time density (FPTD) fi(t), we 
have to consider three possible cases that may occur for 
each non-default firm i (for more details, we refer to ref- 
erences [H] and [T3]): 

1. First passage happens inside the interval. We 

know that if X l (T+_ 1 ) > A(Tj_i) and Xi(Tf) < 
Di(Tj), then the first passage happened in the time 
interval [Tj_i,Tj]. To evaluate when the first pas- 
sage happened, we introduce a new viable bij as 
r ' Tj ~ 1 We generate several correlated uni- 



bij = 



form numbers Yi by using the SOU method, then 
compute Si — bijYi + Tj_\. If s% belongs to interval 
[Tj_i, Tj], then the first passage time occurred in this 
interval. We set IsDef ault(i) = 1 to indicate firm 
i has defaulted. To get the density for the entire 

interval [0,T], we use /,,„(*) = ( ^{Z^' j 9ij(si) * 
K(h opt ,t — Si), where n is the iteration number of 
the Monte-Carlo cycle, and gij(si) is the conditional 
boundary crossing density used to obtain an appro- 
priate density estimate [I]. 

2. First passage does not happen in this inter- 
val. If Si does not belong to interval [Tj_±,Tj], then 
the first passage time has not yet occurred in this 
interval. 

3. First passage happens at the right boundary 
of the interval. If Xi(Tf) < AC 7 }) and X^Tf) > 
Di(Tj), then Tj is the first passage time. We evaluate 
the density function using kernel function /i jU (t) = 
K(h opt ,t — Tj), and set IsDefault(i) = 1. 

Next, we increase j and examine the next interval and 
analyze the above three cases for each non-default firm 
again. After running N times Monte-Carlo cycle, we get 
the FPTD of firm i as %(t) = 4 £^ =1 %, n (t). 



4 Generalizations 

In the above algorithms, we only consider one driving 
Poisson shock that affecting all the firms. Sometimes it is 
necessary to have several independent shocks to account 
for events that affect individual companies rather than 
the entire market. Hence, our next goal is to generalize 
the developed multivariate uniform sampling (MUNIF) 
method to the case of several independent Poisson shocks. 

We consider d firms which are driven by m independent 
Poisson shocks N£ , • • • , N™ as described by Eq. ^ . Let 
Mfc be the number of jumps for each Poisson shock Nf, 



M = Y^k=\ is the total number of jumps. We gen- 
erate the jump instant Tj^ by generating the interjump 
times (Tj if — Tj—%^) for each Poisson shock. Then we sort 
the jump instant Tj k by the relevant ascending order and 
still denote them as Tj (j = 1, 2, • • • , M). Furthermore, 
an array ShockType (whose size is M) is used to record 
the type of shock at Tj . Then we can carry out multivari- 
ate uniform sampling method for this case as same as in 
Section [31 but the postjump value should be calculated 
as: 

X i (T+)=X i (Tf)+Y i (T j ), 

where Yi(Tj) is the size of jump for i-th firm at Tj, the 
type of shock is determined by the array ShockType. 
Besides, we may generate correlated Yi(Tj) for different 
firms. 



5 Applications and discussions 

In this section, we demonstrate the developed model at 
work for analyzing the default events of two differently 
rated and correlated firms (denoted as A- and Ba-rated 
firms according to the Moody's debt rating system) via a 
set of historical default data as presented by [5] . Our first 
task is to describe the first passage time density functions 
and default rates of these firms. 

In order to apply our developed procedure, first we need 
to calibrate the developed model, in other words, to nu- 
merically choose or optimize the parameters, such as drift, 
volatility and jumps to fit the most liquid market data. 
As mentioned in Section [3l after Monte-Carlo simulation 
we obtain the estimated density fi(t) by using the kernel 
estimator method. The cumulative default rates for firm 
i in our model is defined as, 



fi(r)dT. 



(9) 



Then we minimize the difference between our model and 
historical default data Ai(t) to obtain the optimized pa- 
rameters in the model (such as cr^- , arrival intensity A in 
Eq. O): 



argmin 



\ 



E 



Pi(tj)-Mtj) 



(10) 



For convenience, we reduce the number of optimizing pa- 
rameters by: 

1. Setting X(0) = 2 and ln(/t) = 0. 

2. Setting the growth rate 7 of debt value equivalent to 
the growth rate \i of the firm's value [5], so the default 
of firm is non-sensitive to /1. In our computations, 
we set fi = — 0.001 . 

3. The interjump times (Tj —Tj-\) satisfy an exponen- 
tial distribution with mean value equals to 1. 

4. The arrival rate for jumps satisfies the Poisson distri- 
bution with intensity parameter A, where the jump 
size is a normal distribution Z t ~ N(fiz, &z)- 



As a result, we only need to optimize a, A, pz, ®z 
for each firm. This is done by minimizing the difference 
between our simulated default rates and historical data. 
Moreover, we assume there is only one driving Poisson 
shock with arrival rate A. Hence we first optimize four 
parameters for, e.g., the A-rated firm, and then set the 
same A for Ba-rated firm. 

The minimization was performed by the using quasi- 
Newton procedure implemented as a Scilab program. The 
optimized parameters are described in Table [1] 

Table 1: Optimized parameters for A- and Ba-rated firms 
by using the MUNIF method. In each step of the opti- 
mization, we choose the Monte Carlo runs N = 50, 000. 
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Figure 3: Density function (top) and default rate (bot- 
tom) of Ba-rated firm. The simulations were performed 
with Monte-Carlo runs TV = 100, 000, for the conventional 
Monte-Carlo method, the discretization size of time hori- 
zon was A = 0.005. 



By using these optimized parameters, we carried out the 
final simulation with Monte Carlo runs N = 100, 000. 
Moreover, we also carried out simulations based on the 
conventional Monte-Carlo method with the same parame- 
ters and the discretization size of time horizon A = 0.005. 
The estimated first passage time density functions of A- 
and Ba-rated firms are shown in the top of Fig. O and 
[31 respectively. The simulated cumulative default rates 
(line) together with historical data (squares) are given in 
the bottom of Fig. [2] and 03 In Tabled we give the cal- 
culated optimal bandwidths and the corresponding CPU 
times. Based on these results, we conclude that: 




Figure 2: Density function (top) and default rate (bot- 
tom) of A-rated firm. The simulations were performed 
with Monte-Carlo runs N — 100, 000, for the conventional 
Monte-Carlo method, the discretization size of time hori- 
zon was A = 0.005. 

1. A-rated firm has a smaller Brownian motion part 
compared with Ba-rated firm, besides B-rated firm 
has large pz and especially large az, which indicate 
that the loss due to sudden economic hazard may 
fluctuate a lot for Ba-rated firms. 

2. The density function of A-rated firm still has the 
trend to increase, which means the default rate of A- 
rated firm may increase little faster in future. As for 
Ba-rated firm, its density functions has decreased, so 



Table 2: The optimal bandwidth h op t and CPU time per 
Monte-Carlo run of the simulations. The simulations were 
performed with Monte-Carlo runs N = 100, 000, for the 
conventional Monte-Carlo (CMC) method, the discretiza- 
tion size of time horizon was A = 0.005. 
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its default rate may increase very slowly or be kept 
at a constant level. 

3. From the CPU time in Table we can conclude 
that the multivariate UNIF approach is much more 
efficient compared to the conventional Monte-Carlo 
method. 

Our final example concerns with the default correlation 
of the two firms. We use the following conditions in our 
multivariate UNIF method: 

1. Setting X(0) = 2 and La(/s) = for all firms. 

2. Setting j = p and p = —0.001 for all firms. 

3. Since we arc considering two correlated firms, we 
choose a as, 



and 
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(11) 



(12) 



In Eq. (|12p. p reflects the correlation of diffusion 
parts of the state vectors of the two firms, we set 
p = 0.4 as in [5j. Furthermore, we use the optimized 
<7i and (72 in Table [1] for A- and Ba-rated firms, re- 



spectively. Assuming a 12 = 0, we get, 



6 Conclusions 




0, 

P12CT2, 



^22 = V 1 Pl2 CT 2- 



4. The arrival rate for jumps satisfies the Poisson dis- 
tribution with intensity parameter A = 0.1 for both 
firms. The jump size is a normal distribution Z t ~ 
N{nznVZi), where fiZi and cr^ can be different for 
different firms to reflect specifics of the jump process 
for each firm. We adopt the optimized parameters 
given in Table [TJ 

5. As before, we generate the same interjump times 
(Tj — Tj-i) that satisfy an exponential distribution 
with mean value equals to 1. 



In this contribution, we developed efficient Monte-Carlo- 
based computational procedures for the solution of the 
FPT problem in the context of multivariate (and cor- 
related) jump-diffusion processes. This was achieved 
by combining a fast Monte-Carlo method for one- 
dimensional jump-diffusion process and the generation 
of correlated multidimensional variables. The devel- 
oped procedures were applied to the analysis of mul- 
tivariate and correlated jump-diffusion processes. We 
have also discussed the implementation of the devel- 
oped Monte-Carlo-based technique for multivariate jump- 
diffusion processes driving by several compound Poisson 
shocks. Finally, we have applied the developed technique 
to analyze the default events of two correlated firms via 
a set of historical default data. The developed methodol- 
ogy provides an efficient computational technique that is 
applicable in other areas of credit risk and pricing options. 



The simulated default correlations can be obtained via 
the following formula: 
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(13) 



where Pi2,n is the probability of joint default for firms 1 
and 2 in each Monte Carlo cycle, P\. n and P 2 . n are the 
cumulative default rates of firm 1 and 2, respectively, in 
each Monte Carlo cycle. 

The simulated default correlations of A- and Ba-rate 
firms are given in Fig. 0J Based on these results, we can 
conclude that 

1. The default correlations tend to increase over long 
horizons and may converge to a stable value. 

2. Our developed methodology gives almost identical 
default correlations compared with the conventional 
Monte-Carlo method which confirms the validity of 
the developed methodology. 
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Carlo method, the discretization size of time horizon was 
A = 0.005. 
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